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Abstract 



Initial data corresponding to spacetimes containing black holes are considered 
in the time symmetric case. The solutions are obtained by matching across 
the apparent horizon different, conformally flat, spatial metrics. The exterior 
metric is the vacuum solution obtained by the well known conformal imaging 
method. The interior metric for every black hole is regular everywhere and 
corresponds to a positive energy density. The resulting matched solutions 
cover then the whole initial (Cauchy) hypersurface, without any singularity, 
and can be useful for numerical applications. The simpler cases of one black 
hole (Schwarzschild data) or two identical black holes (Misner data) are ex- 
plicitly solved. A procedure for extending this construction to the multiple 
black hole case is also given, and it is shown to work for all time symmetric 
vacuum solutions obtained by the conformal imaging method. The numerical 
evolution of one such 'stuffed' black hole is compared with that of a pure 
vacuum or 'plain' black hole in the spherically symmetric case. 

PACS numbers: 04.70.Bw,04.20.Cv 
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I. INTRODUCTION AND OVERVIEW 



Black holes are the most elementary objects in General Relativity. Allowing for the well 
known 'no hair' theorems, they can be characterized by their mass, charge and spin, like 
elementary particles in classical or quantum mechanics. This is because event horizons act 
as a one way membrane which separates the black hole interior from the exterior: there is 
no causal effect that can propagate from the inner to the outer region. Therefore, one does 
not need to know all the details about physical processes taking place at the interior in order 
to describe the overall evolution of the black hole. This also means that black hole metrics 
that differ only inside the holes can lead to the same exterior spacetime. The idea of this 
paper is to replace the singular vacuum solution for a black hole interior by a regular one 
corresponding to a non-vanishing energy density, but keeping the same exterior metric. The 
term 'stuffed' black hole is a good description of the resulting solution. 

The 'no hair' property has been very useful in Numerical Relativity to avoid the interior 
black hole singularity appearing in the initial data. The standard practice is to excise the 
inner region from the computational domain, so that one can safely evolve the exterior 
(usually vacuum) region. This has been done for a single black hole (Schwarzschild case) by 
setting an internal boundary on the computational domain at the initial position tq of the 
apparent horizon. The same thing has been done for a system consisting of two non rotating 
black holes, starting the evolution from the Misner data or even for an arbitrary number 
of them, by using the 'conformal imaging formalism' |3| to construct conformally flat initial 
data. We shall refer to all these pure vacuum solutions as 'plain' black hole metrics to 
distinguish them from the 'stuffed' ones we are presenting here. 

In this paper we propose to take advantage of the 'no hair' property in a different 
way. We will keep the same exterior metric of the 'plain' black hole case, obtained by 
the conformal imaging formalism. But we will match this exterior metric to a regular 
conformally flat interior corresponding to a positive energy density. The exterior solution is 
not affected because the matching surfaces coincide in every case with the initial position 
tq of the black hole apparent horizon. This is done for a single non rotating black hole 
(Schwarzschild metric) in Section II. The interior metric in this case is the spatial part 
of a positive curvature Friedmann-Robertson- Walker (FRW) metric, so that the resulting 
construction can be interpreted as initial data for the Oppenheimer-Schneider dust collapse 
(Section III). In Section IV, we give in closed form the corresponding solution for two non- 
rotating stuffed black holes (the analogous of the Misner data). We provide in Section V 
that a step by step procedure to stuff time symmetric initial data containing any number of 
black holes with arbitrary sizes and locations. 

Finally, in the Appendix, we compare the numerical evolution of a 'plain' black hole 
with that of a 'stuffed' black hole in the spherically symmetric case. This is a biased test, 
because spherical coordinates are singular at the origin and this coordinate singularity affects 
only the stuffed hole, because in that case there is no excised region and the computational 
domain includes the origin. In spite of that, we have evolved stuffed holes with the same 
accuracy and stability as plain holes. This opens the door to three-dimensional applications, 
where the 'stuffing' approach avoids putting internal curvilinear boundaries that would need 
a special treatment. The advantage of the 'stuffed' black holes versus the 'plain' ones will 
be manifest if the coordinate system contains a shift vector which allows the holes to move 
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across the numerical grid: in the plain case one would need to treat a number of moving 
curvilinear boundaries, but in the stuffed case nothing special is to be done because there 
are no internal boundaries in that case. 

II. SINGLE HOLE INITIAL DATA 

The evolution formalism is specially well suited for Numerical Relativity. In normal 
coordinates, the spacetime metric can be written as 

ds 2 = —a 2 dt 2 + jij dx % dx 1 . (1) 

To start the evolution one needs both the spatial metric 7^ on the initial hypersurface and 
its second fundamental form (extrinsic curvature) 

K a = dt ^i • ( 2 ) 

These initial data, however, are constrained by the following equations 

16tt t = 16tt a 2 T 00 = R J j + (K)) 2 - K\K) (3) 
(energy constraint), and 

8tt Si = 8tt a T° = K j i;j - d^K)) (4) 

(momentum constraint), where both the Ricci tensor and the covariant derivatives are the 
ones corresponding to the three-dimensional geometry of the initial slice. These constraint 
equations are first integrals of the Einstein equations: they are verified at every time slice 
by the spatial part of any spacetime metric. 

One can learn a lot about the spacetime by just looking at the geometry of the spatial 
hyper surf aces. Let us consider for instance a closed two-surface a in one such a three- 
dimensional manifold. The expansion 9 of a congruence of outgoing light rays starting at a 
is given by 

9 = ir j :j + Kij n* n j - K\ (5) 

where the three-dimensional vector n l is the unit normal to a. In the 'time symmetric' case 
(K^ = 0) the minimal surfaces of the three- geometry 

n j ;j = (6) 

are also apparent horizons (9 = 0) and vice versa. 

In the case of non-rotating black holes, one usually starts with time symmetric initial 
data (and therefore zero initial momentum density Si) and conformally flat initial metrics, 
namely 

la = ^ k ■ ( 7 ) 



.3 



This simplifies the constraint equations ||, so that (|3|) can be written in terms of a flat 
space Laplace operator A 



A ijj = -2tt t ip 5 . (8) 

In the spherically symmetric case, the regular vacuum solution of @ can be written as 

in 

?/>sh/ = 1 + - , (9) 

which corresponds to the space part of the Schwarzschild metric in isotropic coordinates. 
Allowing for (|7|), the minimal surface (apparent horizon) condition (P) reads 

d r ( r ^ 2 ) = 0, (10) 

which holds for = m/2. Note that the combination r ip 2 is precisely the 'area radius' R 
of the spherical geometry, so that the apparent horizon (minimal surface) is placed at 

R(r ) = 2m . (11) 

Condition ([H]) can also be obtained by realizing that the metric given by (|7|||) is invariant 
under the discrete symmetry 

r < > Tq/t . (12) 



which can be interpreted as a coordinate inversion at r = tq. This inversion symmetry 
provides suitable boundary conditions when one excises the spherical region r < tq from the 
computational domain. The resulting inner boundary condition is very easy to implement 
numerically in spherical coordinates, but it is much more difficult to manage in the three- 
dimensional case (where one would get a curvilinear boundary in cartesian coordinates) even 
for a single black hole. 

As an alternative to this inner boundary aproach, we will construct a 'stuffed' black hole 
by replacing the initial data (§) with 

1 + for r > m/2 

,p = < Ittw tmr<m ' 2 <13) 

so that in the exterior region one recovers the 'plain' black hole solution (P), whereas in the 
interior region one gets, allowing for (|7|), a homogeneous three-dimensional metric which is 
the space part of a closed Friedman- Robertson- Walker (FRW) model (see Fig. 1). 



III. PHYSICAL INTERPRETATION 

Up to now, we have just constructed time symmetric data on the initial slice. This 
means that we have used only the energy and momentum constraints, but not the remaining 
(evolution) field equations. From the energy constraint (|8|), one can compute the energy 
density 
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{0 for r > m/2 

3 2 for r < m/2 ( 14 ) 
3277 77T ' 

whereas the momentum constraint (U) plus the time-symmetry condition (fQj = 0) imply 
that the momentum density should vanish. However, in order to get a physical interpretation 
of the matter content of the spacetimes generated by these initial data, we need to say 
something about the remaining (space) components of the stress-energy tensor T^ w . 

To do this, we can study the time evolution of the matching conditions between the 
constant density interior and the vacuum exterior, namely 

$^ =0 = , (15) 

where $ = is the equation of the matching hypersurface. Allowing for the fact that the 
momentum density vanishes, one gets easily from fljjD that $ can not depend on time. This 
means that the matching surface can be taken to be the sphere 

$ = r - m/2 (16) 



for any value of time. It follows then from (|15|) that the radial direction gives an eigenvector 
of the stress-energy tensor (computed at the matching surface) with zero eigenvalue. 

If we assume for simplicity a barotropic perfect fluid matter content, it follows that the 
pressure should vanish (incoherent matter). The initial data ( pT3|) can then be understood as 
corresponding to a particular case of the well known Oppenheimer-Schneider dust collapse: 
a constant density spherical star which is initially at rest. In our case, the initial star radius 
actually coincides with the position of the apparent horizon. 

The former is just one of the many possibilities allowed by the no-hair theorems. A very 
interesting case is the 'string perfect dust' matter content || 

= S» p S» , (17) 

where 

is a simple surface-forming bivector. 

There are basically two different possibilities to get a zero radial eigenvalue in the string 
case [[|: either S* Mr = (the bivector has no radial component) or the vector A 1 ' = is 
isotropic (and, being antisymmetric, has no radial component). In both ways, it follows that 
the string stress tensor breaks the spherical symmetry of the initial data, and this will lead 
to a non-spherical time evolution. This 'dynamical symmetry breaking' is in contrast with 
the more predictable behaviour that one gets in the dust case, where the full stress-energy 
tensor is spherically symmetric, so that the spherical symmetry of the initial data will be 
preserved during evolution. 



IV. TWO BLACK HOLES INITIAL DATA 

The Misner initial data [0 are axially symmetric and describe two identical non-rotating 
black holes which are initially at rest. They can be obtained by linear superposition of 
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spherically symmetric solutions of the Laplace equation, with centers distributed along the 
symmetry axis: 



^ = ! + £($+ + $;) (19) 

n=l 

with 

*n = "T csch(n// ) , (20) 



where 



( r n ) 2 = a; 2 + y 2 + [2 ± a coth(n/i )] 2 . (21) 



In order to see how this solution is obtained, let us notice that in the time symmetric 
vacuum case the constraint (|8|) on the conformal factor ip is the flat space Laplace equation. 
We know from electrostatics that we can take advantage of the invariance of the Laplace 



equation under discrete symmetries, such as (|12|), by using the 'imaging method' [10]. This 



technique was adapted to Relativity by Misner ]12| in order to obtain initial data which 



are invariant by inversion across a number of spherical surfaces, which will become at the 
end minimal surfaces (apparent horizons) of the resulting solution. This solution can be 
interpreted as describing time symmetric initial data for a number of black holes. 
The apparent horizons (minimal surfaces) in (|TjD are the two spheres a ± given by 

rf = a csch(/i ) . (22) 

The term has been constructed as the image of $~ under inversion across a + . Also, 
the term has been constructed as the image of under inversion across cr~. This 
means that the dipole combinations 

A± = S± + *ti (23) 

are invariant under inversion across the sphere a^, respectively. 

We have seen in the previous section how to stuff a single black hole. Now we have two 
identical holes and we will look for regular interior solutions for each one. Let us begin by 
considering the interior region to the first sphere a~ . The infinite sum in the Misner solution 
(|19|) can here be written in a more convenient form 

00 

V> = E K (24) 

n=l 

so that every term is now invariant under inversion across the sphere a". We will make 



use of this symmetry property to match every term in fl2jj) to (the conformal factor of) a 
constant curvature metric. 

The conformal factor for a constant (positive) curvature metric (closed FRW) can be 
written as 



if, 



FRW 



(25) 



\ 1 - 2 b 77 cos <p + (b 2 + A 2 /4) (rf ) 



2 > 
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where F is an arbitrary scale factor and A and b are parameters related to the conformal 
transformations of Euclidean three-space. If we impose inversion symmetry (Eq. |T0| ) across 
the sphere a~ , we get 



(6 2 + A 2 /4) a 2 = sinh 2 /i . 



(26) 



This single condition ensures, allowing for fllUp, that the normal derivatives of the vacuum 
and FRW conformal factors will coincide on a~ if and only if both conformal factors actually 
coincide there. 

It follows that, in order to complete the matching, we need only to look at the conformal 
factor values at a~ and tune the free parameters in (|25|) . The FRW conformal factor can be 
easily evaluated on o~: 



FRW\cr- 



F A/2 



\ 1 — b (z — a coth « ) 



where we have noted 

z — a coth /i = t i cos if . 

On the other hand, a straightforward calculation shows that 

2a 



a: 



2$; 



a 2 + 2a z sinh[n/i ] sinh[(n — 1)/xq]/ sinh[/i 



(27) 



(28) 



(29) 



Now it is easy to obtain the values of the arbitrary parameters in (p5|) that ensure the 
matching between the vacuum dipole term A~ and the conformal factor 



F n A r 



\ 1 — 2 b n (z — a coth /i ) + sinh 2 /i (r 1 / a) z 
across cr~ for every value of n. Their actual values are: 

cosh[(2n — l)/io] ~ cosh[/i ] 



(30) 



ab n = — sinh /i 
a X„ 



cosh[(2n — l)/io] cosh[/i ] 
2 sinh[(2n — 1)/U ] sinh 2 [/i ] 



cosh[(2n — l)/io] cosh[/i ] 
4a 



sinh[(2n — l)fio] 



(31) 
(32) 
(33) 



Allowing for the symmetry of the solution across the equatorial plane, it is easy to see 
that the conformal factor 



F n X r 



\ 1 + 2 b n {z + a coth /x ) + sinh 2 /i (r^ / a)' 



(34) 



will also match the vacuum dipole term A+ across a + . It follows that the complete solution 
for the stuffed two hole problem reads 
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J2n=l for r l < a CSch(/i ) 

E£°=i f o r r f < a csch(/i ) (35) 

1 + Y.n=i ( $ n + K ) elsewhere. 



The resulting solution is regular and smooth everywhere (see Fig. 2). Following Bowen 



et al [11], we can define an effective energy density starting from a 'Newtonian potential' 



linearly related to In our case, allowing for (§), this amounts to 

T ef f = T^ 5 . (36) 

This effective energy density can be easily expressed as the sum of the effective energy 
density of every FRW component, namely 



T eff 

where 



En=i( T eff)n for rf < a csch(/i ) ( 37 ) 
elsewhere 



t, 



effjn 



8nF 



«> 5 ■ (38) 



The positivity of the factors ip^ ensures that the energy density fl3~7|) is positive inside the 
holes and it is bounded by the maximum energy density of the FRW components. 

As in the single hole case, the energy density has a jump at the matching surfaces due to 
the discontinuity of the second metric derivatives there. It follows from the previous section 
that a reasonable physical description of the solution could be that of two balls of incoherent 
matter. Notice, however, that the energy density of every ball is not constant, as it was in 
the single hole case. The gravitational interaction between the two balls accounts for their 
distortion. 



V. MULTIPLE BLACK HOLE CASE 

Let us consider now the time symmetric multiple black hole case. As is well known, the 



vacuum exterior solution can be obtained by the conformal imaging method |]12| . As an 
input for this method, one must provide the size and location of an arbitrary number N of 
spheres, which will become at the end the apparent horizons of N black holes. The resulting 
solution, by construction, will be then inversion symmetric across the apparent horizon of 
every black hole. It could be written as a linear superposition of poles, 

* = i + f; a \ . (39) 



n=l 



\x — X 



Our goal is to provide a suitable interior solution for every black hole. As far as the holes 
do not overlap, we can consider them separately. Let us begin with the first one: its horizon 
a* 1 ) is a sphere centered at x%, 

\x-xi\ = Ri. (40) 



S 



Note that, as (^) is invariant under inversion across the sphere cn 1 ), half of the poles in 
are outside and half inside. Moreover, every outside pole is the image under inversion 
of an inside pole and vice versa. 

This means that we can combine every pole with its image to form the invariant dipoles 

"'' ( 41 ) 



+ 



where the exterior point Xj> is the image of the interior one Xj under inversion across the 
sphere The sum in the solution (|39"1) can then be rearranged as follows 



(42) 



Our strategy will be then to match separately every invariant dipole (f4~T|) to a closed 
FRW factor, given by 



(i) 



F 3 X 3 



\l-2b j (x-x 1 ) + {b) + Aj/4) \x-x x \ 



where 



6,- = b 3 



(43) 



(44) 



As in the previous section, we must impose first that the interior solution ( f4~3"D have a 
'throat' fllUP at the sphere obtaining the condition 

(P) + Aj/4) = l/R\ , (45) 

which generalizes (pq). As in the previous section, this condition ensures that the radial 
derivatives of the vacuum and the FRW conformal factors will coincide on the sphere if 
and only if both conformal factors actually coincide there. 
Alowing for fl4"5|), the interior factor (fHf) at the horizon is 



if, 



(i)i 
j la 



Fj Xj/2 



\l- bj (x- 



Xl, 



(46) 



On the other hand, allowing for the fact that the two images in the invariant dipole fl4"T| ) 
coincide by construction on the sphere we have 



A 



(i)i 
j lo- 



rn 



2aj 



R x + 2(x — xi)(xi — xj) + \x\ — Xj 



It follows that the choice of parameters 



bi = -2- 



R l + \xi — Xj\ 



F 3 



2/Ri 


R\- 


fi - 


x j\ 


R\ + 
8 aj 


\x\ — 




Rl + 


\xx - 
9 


Xj\ 





(47) 

(48) 

(49) 
(50) 



ensures the matching between (|4lD and ( ^3|) at 

The interior solution for the first hole can then be obtained by summing the corresponding 
FRW factors (||) 

oo 

*W = 5>f- (51) 
The interior solution for every other hole can be obtained in exactly the same way. 

APPENDIX: ID NUMERICAL EVOLUTION 

We will compare here the numerical evolution of a single spherically symmetric (ID) 
stuffed black hole against a 'plain' one. We use in both cases the same finite difference code 
with a 400 point evenly space numerical grid. The metric quantities and their evolution 
equations are also the same. The matter density and speed are computed in the stuffed 
case by using the standard 'upwind' method to model the continuity equation and the Euler 
equation for dust, respectively (they are equivalent to the stress-energy tensor conservation). 

In the 'plain' case, we have used the throat inversion symmetry to provide the inner 
boundary condition. In the general (3D) stuffed case there is no inner boundary. The only 
concern will arise from the discontinuity of the energy density (second derivatives of the 
metric) at the throat. We have not seen any problem with the matter variables (density and 
speed) in our ID numerical evolution. 

The singularity of spherical coordinates at r = 0, however, demands a special treatment 
of the origin in the ID stuffed case. We have done it by adding two virtual points as the 
mirror image of the first two grid points across the origin. In the 'plain' hole case the origin 
is contained in the excised region and such special treatment is not needed. This is why 
this ID test is biased: it is the best case for the plain hole (inversion symmetry in spherical 
coordinates) and the worst one for the stuffed case (need for a special treatment of the 
origin). 

We have compared the code performance, in terms of accuracy and stability, in both 
cases for three different slicing conditions: harmonic, '1+log' and maximal. In the last two 
cases, we can barely notice any significant difference between the plain and the stuffed hole 
evolution in the exterior region (see Fig. 3). In the harmonic slicing case, however, the code 
for the stuffed black hole crashes as the dust ball in the interior region collapses, evolving 
towards a singularity. This is a consequence of the weaker singularity avoidance properties 
of harmonic slicing. 

These results are promising with regards to 3D applications in rectilinear grids, where 
inversion symmetry across the spherical horizons is not so easy to implement (plain case), 
whereas there is no internal boundary of any kind (the coordinate system is regular every- 
where) in the stuffed case. 
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FIGURES 



FIG. 1. Plot of the conformal factor describing initial data for a single stuffed black hole with 
mass m = 2 (apparent horizon at r = 1 in isotropic coordinates). The smoothness and regularity 
of the solution is evident. 

FIG. 2. Surface plot of the conformal factor describing initial data for two stuffed black holes 
with no = 2, a = 1. The plot is in cylindrical coordinates with the axis along the line joining 
the two holes. The azimuthal angle 4> is suppressed, as the solution is axially symmetric. The 
smoothness and regularity of the solution is evident. 

FIG. 3. Evolution of plain (solid line) and stuffed (dashed line) black holes in the spherically 
symmetric case ('1+log' slicing, 400 grid points). Left: the radial metric coefficient is plotted after 
t = 300m. Right: the maximum error for the black hole mass near the horizon is plotted as a 
function of time. Both the stuffed and the plain case show the same accuracy. 
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